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Abstract. The long time behavior of an initial step resulting in a dispersive 
shock wave (DSW) for the one-dimensional iscntropic Eulcr equations regularized 
by generic, third order dispersion is considered by use of Whitham averaging. 
Under modest assumptions, the jump conditions (DSW locus and speeds) for 
admissible, weak DSWs are characterized and found to depend only upon the 
convexity or concavity of the dispersion and a general pressure law. Two 
mechanisms leading to the breakdown of this simple wave DSW theory for 
sufficiently large jumps are identified: a change in the sign of dispersion, leading 
to gradient catastrophe in the modulation equations, and the loss of genuine 
nonlinearity in the modulation equations. Large amplitude DSWs are constructed 
for several particular dispersive fluids with differing pressure laws modeled by the 
generalized Nonlinear Schrodinger equation. These include supcrffuids (Bose- 
Einstein condensates and ultracold Fcrmions) and "optical fluids". Estimates of 
breaking times for smooth initial data and the long time behavior of the shock 
tube problem are presented. Numerical simulations compare favorably with the 
asymptotic results in the weak to moderate amplitude regimes. Deviations in the 
large amplitude regime are identified with breakdown of the simple wave DSW 
theory. 



PACS numbers: 03.75.Kk, 03.75.Lm, 05.45.Yv, 42.65.Sf, 47.40.Nm 
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1. Introduction 

Nonlinear wave propagation in dispersive media with negligible dissipation can lead to 
the formation of dispersive shock waves (DSWs). In contrast to classical, viscous shock 
waves which are localized, rapid jumps in the fluid's thermodynamic variables, DSWs 
exhibit an expanding oscillatory region bookended by a large amplitude solitary wave 
at one edge and nearly linear waves at the other connecting two disparate fluid states. 
A schematic depicting typical left-going DSWs for positive and negative dispersion 
fluids is shown in figure [TJ These structures are of particular, current interest due 
to their recent observation in superfluidic Bose-Einstein condensates (BECs) of cold 
atomic gases [I1[1I1I1IS] and nonlinear photonics [SllZllHilinillllllllSl. Dispersive 
shock waves also occur in a number of other physical systems including water waves 
[14] (known as undular hydraulic jumps or bores), ion-acoustic plasmas [15] (called 
coUisionless shock waves), and stratified layers in the atmosphere [TU [171 IE! and 
ocean [19]. 

The Whitham averaging technique [101 HI] is the principle analytical tool for the 
dispersive regularization of singularity formation in hyperbolic systems; see e.g., the 
review [22]. For a DSW resulting from step initial conditions, the Riemann problem, 
one can obtain its leading and trailing edge speeds in terms of the left and right 
constant states, analogous to the Rankinc-Hugoniot jump conditions of classical gas 
dynamics. The oscillatory region in space-time inhabited by the DSW is described 
by first order, quasi-linear, modulation partial differential equations (PDEs) for the 
slow variation of the parameters of a nonlinear periodic traveling wave. Whitham 
modulation theory for DSWs was initially developed for intcgrable wave equations 
where the Riemann problem for the hyperbolic modulation equations could be solved 
explicitly for a self-similar, simple wave |231 I24j . An important innovation was 
developed by El [2S] whereby the DSW's trailing and leading edge speeds for the 
Riemann problem could be determined without solving the full set of modulation 
equations. The Whitham-El DSW construction relies on the existence of a simple wave 
solution to the full modulation equations, but does not require its full determination, 
hence analytical results are available even for non-integrable equations. 

In this work, the one-dimensional (ID) isentropic Euler equations arc regularized 
by a class of third order dispersive terms, modeling several of the aforementioned 
physical systems. The time to breaking (gradient catastrophe) for smooth initial 
data is numerically found to fall within bounds predicted by the dispersionlcss Euler 
equations. In order to investigate dynamics post-breaking, the long time resolution of 
the Riemann problem is considered. The DSW locus relating upstream, downstream 
flow configurations and the DSW speeds for admissible weak shocks are determined 
explicitly for generic, third order dispersive perturbations. The results depend only 
upon the sign of the dispersion {sgnojkk) and the general pressure law assumed. A 
fundamental assumption in the DSW construction is the existence of an integral curve 
(simple wave) of the Whitham modulation equations connecting the upstream and 
downstream states in an averaged sense. Explicit, verifiable criteria for the breakdown 
of the simple wave assumption arc given. The regularization for large amplitude 
DSWs depends upon the particular form of the dispersion. Thus, DSWs are explicitly 
constructed for particular pressure laws and dispersive terms of physical origin 
including a generalized Nonlinear Schrodingcr (gNLS) equation modeling supcrfiuids 
and nonlinear optics. Comparisons with a dissipative regularization are presented in 
order to highlight the differences between viscous and dispersive shock waves. 
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s_ = X v_ ^ ^ f + X s + = 



Figure 1. Density for the negative dispersion, 1~-DSW case (a) and positive 
dispersion, 1+-DSW case (b) with stationary soliton edge s+ = s_ = (see 
section [6ll. The background flow velocities ui, U2 and linear wave edge velocities 
v^, V— are also pictured. In (a), backflow < 0) occurs while in (b), it is 
possible for the downstream flow to be negative when a vacuum point appears 
(see section |9ll. 

The outline of this work is as follows. Section [5] presents the general dispersive 
Euler model and assumptions to be considered, followed by section [3] outlining the 
gNLS and other dispersive Eulerian fluid models. Background sections [4] and [5] review 
the theory of the hyperbolic, dispersionless system and the Whitham-El method of 
DSW construction, respectively. A detailed analysis of DSW admissibility criteria and 
the breakdown of the simple wave assumption are undertaken in section [5] followed by 
the complete characterization of admissible weak DSWs in section [T] The theory of 
the breaking time for the gNLS model is shown to agree with numerical computation 
in section [H Large amplitude DSWs for the gNLS equation are studied in section 
[9l The manuscript is completed by conclusions and an appendix on the numerical 
methods utilized. 

2. Dispersive Euler Equations and Assumptions 

The ID dispersive Euler equations considered in this work are, in non-dimensional 
form 

Pt + [pu),, = 0, 

{pu)t + [pu^ + -P(p)]^ = [D{p,u)]x, -oo < X < oo, 

where p is a fluid density, u is the velocity, and is the (conservative) dispersive 
term. Formally setting D = gives the hydrodynamic approximation, valid until 
gradient catastrophe when the dispersion acts to regularize the singular behavior. 
The dispersionless, hyperbolic, isentropic Euler equations are known as the P-system 
whose weak solutions to the Riemann problem are well-known |261 127] . Here, the long 
time behavior of the dispersively regularized Riemann problem is analyzed. By the 
formal rescaling X = sx, T = st, 

Pt + {pu)x = 0, 

{pu)t + [pu^ + P{p)\^ = e^[D{p,u)]x, -oo < X < oo,' 

the long time (t ^ 1) behavior of the dispersive Euler equations in the independent 
variables {X,T) is recast as a small dispersion (e^ <C 1) problem. Due to the 



Shock waves in dispersive Eulerian fluids 4 

oscillatory nature of the small dispersion limit, it is necessarily a weak limit as 
shown rigorously by Lax, Levermore, and Venakides for the Korteweg-deVries equation 
(KdV) [251 Uni 1301 131] • In this work, the multiscale Whitham averaging technique will 
be used to study the behavior of the dispersive Eulcr equations (|2.ip for sufficiently 
large time and long waves. 

The Whitham-El DSW simple wave closure method [25] is used to construct 
DSWs under the following assumptions. 

Al (sound speed) The pressure law P = P{p) is a smooth, monotonically increasing 
function of p, P'{p) > for p > so that the speed of sound 



Co = c{po) = y/P'ipo), 

is real and the local Mach number 

^"0 = —, — r, 
c(Po) 

is well-defined. It will also be assumed that the pressure is convex P"{p) > for 
p > so that c'{p) > 0. 

A2 (symmetries) Ec^uations (|2.ip admit the Galilean invariance: 

D{p, u — uq){x — Uf)t, t) ~ D{p, u){x, t), 
for all Mo G M and exhibit the sign inversion 

D{p, -u){^x, t) = -Dip, u){x, t), (2.2) 
so that (|2.1|) arc invariant with respect to x — > —x, m— > —u. 

A3 (dispersive operator) The dispersive term {D[p,u])x is a differential operator 
with D of second order in spatial and/or mixed partial derivatives such that the 
system p.2p has the real-valued dispersion relation 

w = uofc±a;o(fc,/9o), (2.3) 
with two branches found by linearizing about the uniform background state 
p = pq, u = uq with small amplitude waves proportional to exp[i(fc.T — ^^t)]. 
The appropriate branch of the dispersion relation is fixed by the ± in (j2.3|) with 
ujo{k,po) > ior k > 0, po > 0. The dispersion relation has the long wave 
expansion 

uj„{k,po) = cok + pk^ + o{k^), fc -> 0, 7^ 0. (2.4) 

The sign of the dispersion is sgn Wq (fc; po) for fc > 0. Using (|2.4p and the convexity 
or concavity of ujq as a function of fc, one finds 

sgn(^^^j^=sgn^(fc,po). 

Therefore, positive dispersion corresponds to increasing phase and group velocities 
with increasing k while negative dispersion leads to decreasing phase and group 
velocities. 

A4 (Whitham averaging) Equations ()2.ip are amenable to Whitham averaging 
whereby a DSW can be described by a slowly varying, single-phase traveling 
wave. This requires 

(i) The system possesses at least three conservation laws. The mass and 
momentum equations in (|2.ip account for two. An additional conserved 
quantity is required. 
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(ii) There exists a four parameter family of periodic traveling waves parametrized 
by, for example, the wave amplitude a, the wavenumber k, the average 
density p, and the average velocity u limiting to a trigonometric wave for 
small amplitude and a solitary wave for small wavenumber. In the cases 
considered here, the periodic traveling wave manifests as a solution of the 
ordinary differential equation (ODE) (p')^ = G{p) where G is smooth as it 
varies over three simple, real roots. Two roots coincide in the small amplitude 
and solitary wave limits. 

A5 (Simple wave) The Whitham-El method requires the existence of a self-similar 
simple wave solution to the four Whitham modulation equations (the averaged 
conservation laws and the conservation of waves). For this, the modulation 
equations must be strictly hyperbolic and genuinely nonlinear. 

Assumption Al provides for a modulationally stable, hydrodynamic long wave 
limit. The symmetry assumptions in A2 are for convenience and could be neglected. 
As will be demonstrated in section [3l A3 is a reasonable restriction still allowing for a 
number of physically relevant dispersive fluid models. The assumptions in A4 and A5 
allow for the application of the Whitham-El method. While the assumptions in A4 are 
usually verifiable, A5 is often assumed. Causes of the breakdown of assumptions A3 
(unique dispersion sign) and A5 (genuine nonlinearity) are identified and associated 
with extrema in the DSW speeds as either the left or right density is varied. 

The nonstationary DSW considered here is the long time resolution of an initial 
jump in the fluid density and velocity, the Riemann problem 

/ \ ui X <Q , . I pi a; < , , 

u{x,0) = < , p{x,0) = < , (2.5) 

\^ U2 X > P2 X > 

where Uj £ R, pj > 0. 



3. Example Dispersive Fluids 

The dispersive Euler equations (j2.2p model a number of dispersive fluids including, 
among others, superfluids and optical fluids. The particular model equations described 
below were chosen because they incorporate different pressure laws and allow for 
different signs of the dispersion, key distinguishing features of Eulerian dispersive 
fluids and their weak dispersion regularization. 



3.1. gNLS Equation 

The generalized, defocusing nonlinear Schrodinger equation 

/(0)=0, /(p)>0, p>0, 

or gNLS, describes a number of physical systems. For example, the "polytropie 
superfluid" 

f{p) = pP, p> 0, (3.2) 

corresponds to the cubic NLS when p — 1 that describes a repulsive BEC and intense 
laser propagation through optically defocusing (normal dispersion) media. The model 



Shock waves in dispersive Eulerian fluids 



6 



(|3.2[) with p = 2/3 describes a zero temperature Fermi gas near unitarity [311 133] 
which is of special significance as recent experiments have been successfully interpreted 
with both dissipative [31| and dispersive [35] regularizations. Moreover, the regime 
2/3 < p < 1 describes the so-called BEC-BCS transition in ultracold Fermi gases |36j . 
The quintic NLS case, p — 2, models three-body interactions in a BEC [ST] [38]. A BEC 
confined to a cigar shaped trap exhibits effective ID behavior that is well-described 
by the non-polynomial nonlinearity |391 140| 



f{p) = , 7 > 0, (3.3) 

7 

here scaled so that f{p) — >■ p, 7 — >■ 0+. In spatial nonlinear optics, photorefractive 
media corresponding to [HJ |35] 

f{p) = 7 > 0, (3.4) 

1 + 7P 

is of particular interest due to recent experiments exhibiting DS Ws [S] [HI 13 [131 [M] ■ 
For < 7 <C 1, the leading order behavior of p.3p and p.4p correspond to the cubic 
NLS. 

The complex wavefunction -0 can be interpreted in the dispersive fluid context by 
use of the Madelung transformation [45] 

V' = ^/pc'^ (3.5) 

Using (j3.5p in p.ip and equating real and imaginary parts results in the dispersive 
Euler equations (|2.1|) with 

p 



Pip) = / pfmp, c{p) = vvno), 

(3.6) 



[i?(p,u)]. = i[p(lnp)J^ = ^ 



The dispersive regularization of ()2.ip corresponds to the semi-classical limit of p.l|) . 
which, in dimensional units corresponds to ?i — >■ for quantum many body systems. In 
applications, the dispersive regularization coincides with a strongly interacting BEC 
or a large input optical intensity. 

Assumption Al restricts the admissible nonlinearity / to those satisfying 

/'(p)>0, (p/'(p))'>0, p>0, (3.7) 

which is realized by p.2p . (|3.3|) generally and for (|3.4p when 7/3 < 1. Assumptions in 
A2 are well-known properties of the gNLS equation [46] . Assumption A3 is clear from 
(|3.6|) and the dispersion relation is 

1 



cjo(fc, p) = fc + fc2/4 - cfc + — |/cKc. (3.8) 

The dispersion is positive because wq^^ (fc; p) > for fc > 0, p > 0. 
Inserting the traveling wave ansatz 

p ~ p(x — Vt), u ~ u{x — Vt), 

into (|2.ip with (|3.6I) and integrating twice leads to 

A 

u^V + -, (3.9a) 
P 

42" 



< r f{p)dp+Bp^ 



G{p). {3.9b) 
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It is assumed that G has three real roots pi < P2 5: P3 related to the integration 
constants A, B, and C so that, according to a phase plane analysis, a periodic wave 
exists with maximum and minimum densities p2 and pi, respectively. The fourth 
arbitrary constant, due to Galilean invariance, is the wave speed V. In addition to 
mass and momentum conservation, an additional energy conservation law exists |47| 
which reads 







£t + {u[£ + P(p)]}, = 1 



_ {pu)xPx 

^Pxx 



hence the assumptions in A4 are satisfied. The hyperbolicity of the Whitham equations 
can only be determined by their direct study. The genuine nonlinearity of the 
system will be discussed in section |6l It will be helpful to note the solitary wave 
amplitude/speed relation which results from the boundary conditions for a depression 
(dark) solitary wave 

uo= lim w(^), po= lim p{^), p,„in ee minp(^). 

^|— >cxD >-oo 

A phase plane analysis of (j3.9b|) implies that the roots of G satisfy pi — pmin, 
P2 = Pa = Po resulting in the solitary wave speed s = V satisfying 



(s - Mo) 



po 

{PO - Pmin)/(Po) - / ./(P)dp 



(3.10) 



[PO - Pmin)^ 

The soliton profile can be determined by integration of (j3.9&p . 

Dispersive shock waves for the gNLS equation have been studied for the pure 
NLS case [JHl SS] as well as in ID photorefractive media [5D] and the cubic-quintic 
case [5TJ [S2] . A general DSW analysis will be presented in section [SI 

3.2. Other Systems 

The gNLS equation exhibits positive dispersion. Two additional examples are briefly 
given here with negative dispersion. 

Two-temperature collisionless plasma: The dynamics of the ionic 
component of a two-temperature unmagnetized plasma [S3] satisfy the dispersive Euler 
equations with 

P{p) = p, c{p) = 1, 

D{p,u) ^ -4>l- (pxx, -(j)xx = p - e'^ . 

The electronic potential introduces nonlocal dispersion with dispersion relation 

k k^ 
ujoik,p) ^ ^ -k^—, |fc|<l. 

v/1 + k^/p 2po 

It can be shown that ujkk < 0, fc > thus the system exhibits negative dispersion. 

This system has been analyzed in detail [48] and satisfies assumptions A1-A4. 
Large amplitude dispersive shock waves were constructed in |25| under the assumptions 
of A5. 
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Fully nonlinear shallow water: Shallow waves in an ideal fluid with no 
restriction on amplitude satisfy the generalized Serre equations (also referred to as 
the Su-Gardner or Green-Naghdi equations) [531 IMl [Ml [SZ] with 

p{p) = c{p) = 

1 / 1 A ^^-^^^ 

D{p, ^) = 3 [P^ {"^tx + uuxx -uI)]+(tI pPxx - 2^1) ■ 

The density p corresponds to the free surface fluid height and u is the vertically 
averaged horizontal fluid velocity. The Bond number cr > is proportional to the 
coefficient of surface tension. The dispersion relation is 

1/2 



The sign of the dispersion changes when ujkk = corresponding to the critical values 

p2 ^ , ,1/2 

(Tcr = y or fccr = ~ + 3^1 + P^/c^j 

The critical value (Tcr expresses the fact that shallow water waves with weak surface 
tension effects, a < o'er, exhibit negative dispersion for sufficiently long waves 
(fc < feci) and support elevation solitary wave solutions. Strong surface tension, 
(7 > (Ten corresponds to positive dispersion and can yield depression solitary waves. 
Assumptions A1-A4 hold [57]. DSWs in the case of zero surface tension cr = were 
studied in |58| . 

The properties of DSWs for these systems will be discussed in section [TOl 



4. Background: Dispersionless Limit 



The analysis of DSWs for (|2.ip requires an understanding of the dispersionless limit 

Pt + {pu)x = 0, 

(Ht + K + ^(p)], = o, 



corresponding to f = 0. Equations (|4.ip are the equations of compressible, isentropic 
gas dynamics with pressure law P{p) |59| . They are hyperbolic and diagonalized by 
the Riemann invariants (see e.g. [60] ) 

n^u-f'-^dp', .2=. + /^d,', (4.2) 

with the characteristic velocities 

Xi=u-c{p), \2=u + c{p), (4.3) 



so that 

drj , dr 
ITt 

By monotonicity of 



^ + Aj^ = 0. , = 1,2. ,4.4) 



Shock waves in dispersive Eulerian fluids 



the inversion of 



is achieved via 



1 



" = o(''l + ^2), 



(7'2 - ri) 



(4.5) 



2' w, r 

In what follows, an overview of the properties of equations (j4.1[) is provided for 
both the required analysis of DSWs and for the comparison of classical and dispersive 
shock waves. 



4-1. Breaking Time 

Smooth initial data may develop a singularity in finite time. The existence of Riemann 
invariants (|4.2p allows for estimates of the breaking time at which this occurs. In what 
follows, Lax's breaking time estimates [SI] are applied to the system (|4.4p with smooth 
initial data. 

Lax's general approach for 2x2 hyperbolic systems is to reduce the Riemann 
invariant system (|4.4p to the equation 

z' = -a(t)z2, z(0) = m, (4.6) 

along a characteristic family, ' = ^ + A^^, and then bound the breaking time by 
comparison with an autonomous equation via estimates for a and m in terms of initial 
data for ri, r2- 

Following Lax [61], integration along the 1-characteristic family in (14. 6p leads to 
z = e'^dri/dx, a = e~'^dXi/dri and h{ri,r2) satisfies 



dh 
dr2 



9Ai 
dr2 



Ai — A2 

By direct computation with eqs. (|4.2I) . (|4.3p . one can verify the following 



1 



■In 



c(p) 



_^d\i _ c{p)+pc'{p) 



dri 
dx 



dri 
dx 



Hp) 

'c{p) 



cip) 



1/2 



1/2 



(4.7a) 
(4.7&) 
(4.7c) 



The initial data for ri and r2 arc assumed to be smooth and bounded so that they 
satisfy 

ri<ri{x,t) r2<r2(x,t)<r^, 0<t<t^,. (4.8) 

Assuming p > (non- vacuum conditions) , then (|4.76P implies a > and z is decreasing 
along the 1-characteristic. Bounds for a{t) are defined as follows 

B = max a, (4.9) 

peRE ■ 



A 



mm a, 

peRA 



where Ra and Rb are intervals related to the bounds on the initial data, chosen 
shortly. The initial condition m is chosen as negative as possible 

du c{p) dp'\ [c(p)"'^^^ 
dx 



xn = argmin2(0) = argmin 



m = min z(0) 



du c{p) dp 
dx p dx 



p dx 

c{py 



p 



{x,t) = {xo,0) 



(4.10ffl) 
(4.10&) 
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These estimates lead to the following bounds on the breaking time tbr 

- ^ < tbr < ^. (4.11) 

mB mA 

It is still necessary to provide the intervals Ra and Rb- The possible values of 
ri and r2 in (|4.8p and the monotonicity of the transformation for p in (|4.5I) suggests 
taking the full range of possible values Ra = Rb = [9~^{{jj_—rT)/2),g^^{{r2 — ri)/2)]. 
However, this choice does not provide the sharpest estimates in (|4.11l) . The idea is to 
use the fact that ri is constant along 1-characteristics. The choice for m in (j4.10bp 
suggests taking ri = ri(a;o, 0) and allowing r2 to vary across its range of values. While 
the optimal m is associated with this characteristic, it does not necessarily provide 
the optimal estimates for A or i?. A calculation shows 

It can be verified that da/dri < for the example dispersive fluids considered here. 
In this case, any characteristic with ri < ri(a;oiO) can cause a to increase, leading 
to a larger A and a tighter bound in (j4.1ip . If ri > ri(xo,0), then a may decrease, 
leading to a smaller B and a tighter bound on ibr- Combining these deductions leads 
to the choices 

'r2 - ri{xo,0)\ _ifr2-ri'' 



Ra 
Rb 



9 



2 J \ 2 



(4.12a) 
(4.12&) 



when da/dri < 0. 

In summary, given initial data satisfying (|4.8p , the point Xq and m are determined 
from (|4.10aP and (j4.10b| . If m > 0, then there is no breaking. Otherwise, after 
verifying da/dri < 0, the sets Ra and -Rb are defined via (|4.12a|) and ()4.126|) leading 
to A and B in ()4.9p . The breaking time bounds are given by (|4.11|) . A similar 
argument integrating along the 2-characteristic field yields another estimate for the 
breaking time ibr- The only changes are in (|4.10ap and (j4.106p where the minus sign 
goes to a plus sign and the choices for Ra and Rb reflect r2{xo, 0). These results will 
be used to estimate breaking times for dispersive fluids in section El 

4-2. Viscous Shock Waves 

It will be interesting to contrast the behavior of dispersive shock waves for (|2.ip with 
that of classical, viscous shock waves resulting from a dissipative regularization of the 
dispersionless equations. For this, the jump and entropy conditions for shocks are 
summarized below [60l [27l . 

The Riemann problem p.Sp for (|4.ip results in the Hugoniot locus of classical 
shock solutions 

U2=u,±\l^^M^fMiEl^y^\ (4.13) 

I PlP2 J 

The — (+) corresponds to an admissible 1-shock (2-shock) satisfying the Lax entropy 
conditions when the characteristic velocity Ai (A2) decreases across the shock so that 
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P2 > Pi > (pi > P2 > 0). Weak 1-shocks connecting the densities pi and p2 — Pi+A, 
< A ^ 1 exhibit the shock speed 

(i)~wi-ci-i ( — +c; ) A, 0<A<1. (4.14) 



Pi 



Weak, steady (non-propagating) 1-shocks satisfy the jump conditions 

A^^^(Afi-l), 

ci+Pic'i (4.15) 

Ma - 1 - {Ml - 1), < i\fi - 1< 1, 

where Mj = \uj\/cj are the Mach numbers of the up/downstream flows and c'l = 
^(pi). The upstream flow indexed by 1 is supersonic and the downstream flow is 
subsonic, this behavior also holding for arbitrary amplitude shocks. 

4-3. Rarefaction Waves 

Centered rarefaction wave solutions of (|4.1|) exhibit the following wave curves 
connecting the left and right states 

c(p) 

1 — rarefaction : ui = U2 + dp, pi > p2, (4.16a) 

Jpi P 

c{p) 

2 — rarefaction : ui ~ U2 — dp, p-? > pi, (4.16&) 

Jpi P 

where admissibility is opposite to the shock wave case. The characteristic velocities 
Xj increase across a rarefaction wave. Since rarefaction waves are continuous and 
do not involve breaking, the leading order behavior of dispersive and dissipative 
regularizations for (|4.ip are the same. A dispersive regularization of KdV |62j shows 
the development of small amplitude oscillations for the first order singularities at either 
the left or right edge of the rarefaction wave with one decaying as 0{t~^^^) and the 
other C'(t~^/'^). The width of these oscillations expands as 0{t^^^) [23] so that their 
extent vanishes relative to the rarefaction wave expansion with 0(t). 

4.4- Shock Tube Problem 

The general solution of the Riemann problem consists of two waves, each either a 
rarefaction or shock connected by an intermediate, constant state [68l [27] . The shock 
tube problem [59] involves a jump in density for a quiescent fluid ui ~ U2 = 0. The 
solution consists of a shock and rarefaction connected by a constant, intermediate 
state {pm,Um)- For the case pi < p2, a 1-shock connects to a 2-rarefaction via the 
Hugoniot locus (|4.13p (with — ) and the wave curve (|4.16&[) . respectively. For example, 
a polytropic gas with P{p) — np^ gives the two equations 

1/2 



1 — shock 



2 — rarefaction : 



{l^Pla - «Pl)(Pm - Pi) 



PmPl 

2(^7)1/2 



(4.17) 



(7-l)/2 _ (7-l)/2 
P2 Pm 



— 111 O 

7-3 

Equating these two expressions provides an equation for the intermediate density p„ 
and then the intermediate velocity Mm follows. 
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5. Background: Simple DSWs 

The long time behavior of a DSW for the dispersive Euler model ()2.1|) was first 
considered in [25]. In this section, the general Whitham-El construction of a simple 
wave led DSW for step initial data is reviewed. This introduces necessary notation 
and background that will be used in the latter sections of this work. 

Analogous to the terminology for classical shocks, a 1-DSW is associated with the 
Ai characteristic family of the dispcrsionless system (j4.1[) involving left-going waves. 
In this case, the DSW leading edge is defined to be the leftmost (most negative) edge 
whereas the DSW trailing edge is the rightmost edge, these roles being reversed for 
the 2-DSW associated with the A2 characteristic family. The DSW construction for 
1-DSWs is outlined below. A similar procedure holds for 2-DSWs. 

Assuming the existence of a DSW oscillatory region described by slow modulations 
of the periodic traveling wave from assumption A4, three independent conservation 
laws are averaged with the periodic wave. The wave's parameters p, the average 
density, u, the average velocity, fc, the generalized (nonlinear) wavenumber, and a, 
the wave amplitude are assumed to vary slowly in space and time. The averaging 
procedure produces three first order, quasilinear PDEs. This set combined with the 
conservation of waves, fcj + = (w here is the generalized, nonlinear frequency), 
following from consistency of wave modulations, results in a closed system for the 
modulation parameters, the Whitham modulation equations. As originally formulated 
by Gurevich and Pitaevskii [53], the DSW moving boundary value problem is to 
solve the dispcrsionless equations (|4.ip outside the oscillatory region and match this 
behavior to the averaged variables p and u from the Whitham equations at the 
interfaces with the oscillatory region where fc (soliton edge) or a -> (linear 
wave edge). These GP matching conditions correspond to the coalescence of two 
characteristics of the Whitham system at each edge of the DSW. Assumption A5 can 
be used to construct a self-similar, simple wave solution of the modulation equations 
connecting the fc — >■ soliton edge with the a — >■ linear wave edge via an integral 
curve so that the two DSW boundaries asymptotically move with constant speed, the 
speeds of the double characteristics at each edge. In the Whitham-El method, the 
speeds are determined by the following key mathematical observations [25] 

• The four Whitham equations admit exact reductions to quasi-linear systems of 
three equations in the fc — > (soliton edge) and a — > (linear wave edge) regimes. 

• Assuming a simple wave solution of the full Whitham equations, one can integrate 
across the DSW with explicit knowledge only of the reduced systems in the a = 
or fc = planes of parameters, thereby obtaining the DSW leading and trailing 
edge speeds. 

This DSW closure method is appealing because it bypasses the difficult determination 
and solution of the full Whitham equations. Furthermore, it applies to a large class 
of nonintegrable nonlinear wave equations. Some nonintegrable equations studied 
with this method include dispersive Euler equations like ion-acoustic plasmas [53], the 
Serre equations with zero surface tension |58j , the gNLS equation with photorefractive 
|50] and cubic/quintic nonlinearity [63], and other equations including the Miyata- 
Camassa-Choi equations of two-layer fluids [64] . 

Simple DSWs are described by a simple wave solution of the Whitham modulation 
system which necessitates the constancy of one of the Riemann invariants ()4.2[) 
evaluated at the left and right states. Then a necessary condition for a simple DSW 
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is one of 

c(o) 

1- DSW: U2^ui- —dp, P2>Pu (5.1a) 

Jpi P 
f"" dp) 

2- DSW: U2=ui+ —dp, pi > P2- (5.1&) 

Jpi P 

1-DSWs (2-DSWs) arc associated with constant r2 (ri) hence vary in the Ai (A2) 
characteristic field. Equations (|5.1np . (j5.1b[) can be termed DSW loci as they are the 
dispersive shock analogues of the Hugoniot loci (|4.13|) for classical shock waves. It is 
worth pointing out that the DSW loci correspond precisely to the rarefaction wave 
curves in (|4.16a[) and (j4.166p . However, the admissibility criteria for DSWs correspond 
to inadmissible, compressive rarefaction waves where the dispersionless characteristic 
speed decreases across the DSW. The coincidence of rarefaction and shock curves does 
occur in classical hyperbolic systems but is restricted to a specific class, the so-called 
Temple systems [65j to which the dispersionless Euler equations do not belong. 



5.1. Linear Wave Edge 

The integral curve of the Whitham equations in the a = (linear wave edge) plane of 
parameters reduces to the relationships k = k{'p), u = u(j}) and the ODE 

dp u{p) - c(p) - ujk 
where the average velocity is constrained by the density through a generalization of 

u{p) ^ m r ^dp'. (5.3) 



J Pi p 

The negative branch of the linear dispersion relation in p.3p . uj = u(j))k — ujo{k,p), is 
associated with left-going waves, hence a 1-DSW. Using ()5.3|) and ()2.3|) . (|5.2|) simplifies 
to 

dfc ^ cfc/p + a;o_ 
dp c-wo, ■ 

Equation (j5.4p assumes a = 0, an exact reduction of the Whitham equations only at 
the linear wave edge. Global information associated with the simple wave solution 
of the full Whitham equations is obtained from the GP matching condition at the 
soliton edge by noting that the modulation variables satisfy (fc, p, u, a) = (0, pj, Uj, Uj) 
for some j E {1,2} depending on whether the soliton edge is leading or trailing, 
independent of the soliton amplitude aj. Thus (j5.4p can be integrated in the a = 
plane with the initial condition k{pj) = 0, u{pj) = Uj to k{p^-j) associated with the 
linear wave edge, giving the wavenumber of the linear wave edge oscillations. This 
wavepacket's speed is then determined from the group velocity ujk- 

Based on the Riemann data p.5|) . the integration domain for (|5.4|) is p G [pi, P2]. 
The initial condition occurs at either the leading edge where p = pi or the trailing 
edge where p = p2. The determination of the location of the linear wave edge, leading 
or trailing, is set by appropriate admissibility conditions discussed in section |6l The 
solution of ()5.4p with initial condition at pj will be denoted kiji; pj) so that one of 

fc(p,;p,) = 0, .7 = 1,2, (5.5) 
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holds for (j5.4p . Evaluating the solution of (j5.4p at the linear wave edge k{pj; p^^j), 
gives the wavenumber of the linear wavepacket at the leading (trailing) edge when 
J = 1 (j = 2). The associated 1-DSW linear wave edge speed is denoted Wj(pi,P2) 
and is found from the group velocity evaluated at k{pj; ps-j) 

Vj{pi,P2) = ujk[k{pj;p3-j),pj] 

= uiPj) - p3„j), pj], j ==1,2. 

5.2. Soliton Edge 

The speed of the soliton edge is determined in an analogous way by introducing a 
conjugate wavenumber k and frequency 

^{k,-p) = -iw(ifc,p) = u{p)k + \i^o{\k,-p) 

= u{p)k + u}o{k,p). 

The conjugate wavenumber plays the role analogous to an amplitude so that fc — )■ 
corresponds to the linear wave edge where a — > 0. Integrating the ODEs for a simple 
wave in the fc = plane results in 

dfc c(p)fc/p + wo^ 
dp c{p) - 

the same equation as (|5.2p but with conjugate variables. It is remarkable that the 
description of the soliton edge so closely parallels that of the linear wave edge. The 
initial condition is given at the linear wave edge where fc = 0. As in (|5.5|) . the solution 
with zero initial condition at p = pj is denoted fc(p; pj) according to 

fc(pj;p,) = 0, j = l,2. (5.9) 

Then the soliton speed Sj, j = 1,2 is the conjugate phase velocity evaluated at the 
conjugate wavenumber associated with the soliton edge 

,^.(p,,^,)^^fef5zz)lM 

'■(^-^^^^ (5.10) 



k{pf,p3^j) 

Remark 1 In a number of example dispersive fluids studied in this work and 
elsewhere, the transformation to the scaled phase speed 

f-. ujo[k{p),p] 

of the dependent variable in ( [5.^p and the analogous transformation 

c(p)k[p) 

for \5. 8\) are helpful, reducing the ODEs ( |5.^[ ) and \5. 8\) to simpler and, often, 
separable equations for a and a. 
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Figure 2. Integral curves and DSW loci for the dispersive Euler equations with 
c(p) = VP- 

5.3. Dispersive Riemann Problem 

The integral wave curves (|4.16ap . (|4.16b|) and the DSW loci (|5.1ap . (|5.1bp can be 
used to solve the dispersive Riemann problem (|2.5p just as the wave curves and the 
Hugoniot loci are used to solve the classical Riemann problem [27] ■ In both cases, 
the solution consists of two waves, one for each characteristic family, connected by a 
constant intermediate state. Each wave is either a rarefaction or a shock. 

In contrast to the classical case, the integral wave curves and the DSW loci are 
the same for the dispersive Riemann problem. It is the direction in which they are 
traversed that determines admissibility of a rarefaction or a DSW. This enables a 
convenient, graphical description of solutions to the dispersive Riemann problem as 

shown in figure [21 Solid curves ( ) correspond to example 1-wave curves (|5.1ap . 

(|4.16a|) and the dashed curves ( ) correspond to example 2-wave curves (|5.1&|) . 

(j4.16b|) . The arrows provide the direction of increasing dispersionless characteristic 
speed for each wave family. Tracing an integral curve in the direction of increasing 
characteristic speed corresponds to an admissible rarefaction wave. The decreasing 
characteristic speed direction corresponds to an admissible DSW. The solution to 
a dispersive Riemann problem is depicted graphically by tracing appropriate integral 
curves to connect the left state (pi, ui) with the right state (^2,^2). There are multiple 
paths connecting these two states but only one involves two admissible waves. This is 
shown by the thick curve in figure [2] Since the 1-wave curve is traced in the negative 
direction to the intermediate, constant state (pm,Mm), this describes a 1-DSW. The 
2-wave curve is then traced in the positive direction to the right state, describing a 
2-rarefaction. The 1-DSW is admissible because pm > Pi- Since the characteristic 
speed A2 is monotonically increasing, the 2-rarefaction is admissible (pm < 92)- The 
direction of the curve connecting the left and right states was taken from left to 
right. The opposite direction describes an inadmissible 1-rarefaction connected to an 
inadmissible 2-DSW. 

The example shown in figure [5] corresponds to the dispersive shock tube problem 
consisting of an arbitrary jump in density for a quiescent fluid mi = U2 = 0. Such 
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problems have been studied in a number of dispersive fluids, e.g. |66[ 1581 [571 164j . The 
determination of (pmi^m) proceeds by requiring that the left state (pi,0) lie on the 
1-DSW locus ff^THUfi 

um = - r —dp- (5-11) 

Jpl p 

The second wave connects (pm, "i^m) to the right state (/92, 0) via the 2-rarefaction wave 
curve (|4.16&P 

c{p) 



P 



-dp. 



(5.12) 



Equating (|5.1ip and (|5.12p leads to 



1 P 



dp 



P 



dp = 0, 



which determines the intermediate density pi < p^ < p2. The intermediate velocity 
Um < follows. 

This construction of the wave types and the intermediate state (/Omjitni) is 
independent of the sign of dispersion and the details of the dispersive term D in 
(|2.ip . depending only upon the pressure law P(p). For example, polytropic dispersive 
fluids with P{p) = K,p^ (e.g., gNLS with power law nonlinearity p.2p and gSerre 
p. lip ) yield the intermediate state 

2/(7-1) 



j7-l)/2 
Pi 



j7-l)/2 
P2 



1/2 



(5.13) 



3-7 



\,(7-l)/2 ^(7-l)/2 
Pi P2 



This prediction will be compared with numerical computations of gNLS in section 



6. DSW Admissibility Criteria 

As shown in section 14.31 when (|5.1ap holds and pi > P2, sl continuous 1-rarefaction 
wave solution to the dispersionless equations exists. Gradient catastrophe does not 
occur so the rarefaction wave correctly captures the leading order behavior of the 
dispersive regularization. However, when pi < p2, the rarefaction wave solution is 
no longer admissible and integrating the dispersionless equations via the method of 
characteristics results in a multivalued solution. A dispersive regularization leading 
to a DSW is necessitated. Specific criteria are now provided to identify admissible 
1-DSWs. 

The general admissibility criteria for DSWs depend on the ordering of the soliton 
and linear wave edges. In what is termed here a "1"''-DSW", the conditions are HHj 

U2 - C2 < S2 < U2+ C2, (6.1a) 

1+-DSW: vi<ui~ci, (6.1 &) 

S2 > vi, (6.1c) 
where the soliton is at the trailing edge of the DSW. Similarly a "1^-DSW" satisfies 
the conditions 

U2 - C2 < V2 < U2 + C2, (6.2a) 

1"-DSW: si<ui-ci, (6.26) 
V2 > Si, (6.2c) 
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with the sohton at the leading edge. The designation 1+-DSW (1~-DSW) corresponds 
to a positive (negative) dispersion fluid as shown below. These eonditions are 
analogous to the Lax entropy eonditions for dissipatively regularized hyperbolie 
systems [55]. A key difference with classical fluids is that there is only one "sign" 
of dissipation due to time irreversibility. For time-reversible, dispersive fluids, both 
signs are possible. The admissibility criteria ensure that only three characteristics 
impinge upon the three parameter DSW, transferring initial/boundary data into the 
DSW and allowing for the simple wave condition (|5.1a|) to hold. Sufficient conditions 
for these criteria as applied to the dispersive Euler equations are now shown. 

First, consider the criteria for a 1"'"-DSW in a positive dispersion fluid. Inserting 
the linear wave and soliton edge speeds (|5.6p . (|5.10p into inequalities (|6.1a|) . (j6.1bp 
simplifies the first two admissibility criteria to 

- C2 < = < C2, (6.3a) 

ujo^{ki,pi) > ci, (6.3&) 

where k2 = k{p2', pi) and ki = fc(pi;p2)- In section [7.31 the admissibility of weak 
l+-DSWs when < p2 ~ Pi ^ 1 is demonstrated. The further assumptions 

1+-DSW: '^fefe < 0, cfc + pwo^>0, cfc + pwo->0, (6.4) 

enable the extrapolation of Eulerian 1+-DSW admissibility to moderate and large 
jumps p2 > pi (below). Assumptions (|6.4p hold for gNLS, gSerre, and ion-acoustic 
plasma in certain parameter regimes (moderate jumps) . 

The extrapolation of admissibility to larger jumps can be demonstrated as follows. 
For 1-DSWs, the negative branch of the dispersion relation (|2.3p has been chosen for 
A; > 0. Using the small k asymptotics (|2.4p in (|5.4p with initial condition k{p2; P2) = 0, 
one can show that fc(p;p2) is a decreasing function of p for |p — /02I <C P2- Since 
a;Ofc(0,p) = c(p). the convexity of uq implies LijQ^{k,'p) > c{p) for fc > 0. This fact 
combined with (|6.4p in ()5.4p implies that fci = k{pi; P2) > for p2 > pi and that 
inequality ()6.36[) holds. Similar arguments demonstrate that k2 = k{p2\ pi) > for 
P2 > pi and that the inequalities in (|6.3ap hold. Thus, if (|6.4[) hold for an interval 
k G [0,fc*), then (|6.3a[) and ()6.36p are verified for p2 G (pi,/?*) where k{pi; p^) = k^,. 
It is now clear why p2 > pi in (|5.1aP and the designation 1+-DSW is used when the 
sign of dispersion is positive. 

By similar arguments, inequalities (|6.2qP and (j6.2bp hold for negative dispersion 
fluids when pi < p2 and 

1~-DSW : w^j, > 0, ck + puj(^ > 0, ck + pCj^- > 0. (6.5) 

The only change with respect to (|6.4p is the convexity of the conjugate dispersion 
relation. 

The final inequalities ()6.1cP and (|6.2cp require verification. An explicit analysis for 
weak DSWs is given in section [71 An intuitive argument can be given for the general 
case. When uJkk > 0, the ease of positive dispersion, the group velocity of waves 
with shorter wavelengths is larger while the opposite is true for negative dispersion. 
The soliton edge corresponds to the longest wavelength (infinite) hence should be the 
trailing (leading) edge in positive (negative) dispersion systems. For the 1-DSW, a 
trailing soliton edge means s > u as is the case in (|6.1cp . A leading soliton edge 
corresponds to the ordering in (j6.2cp . 

In summary, sufficient conditions for a simple wave led 1-DSW are pi < p2 and 
either (|6.1cp . (|6.4p for a positive dispersion 1+-DSW or (|6.2cp . (|6.5p for a negative 
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dispersion 1 -DSW. Because of this, it is convenient to dispense with the subscripts 
defining the DSW speeds vj and Sj in (|5.6p . (|5.10p and use the notation 

V^=V2, V+=Vi, S- = Si, S+ = S2, (6.6) 

which identifies the dispersion sign. The conditions (|6.1cp or (j6.2cp can be verified a- 
priori while the speed orderings (j6.4[) or (|6.5p must be verified by computing the speeds 
directly. For the case of 2*-DSWs, the requirement is pi > p2 and the appropriate 
ordering of the soliton and linear wave edges. The Lax entropy conditions for the 
dissipativc regularization of the Euler equations give similar criteria, namely positive 
(negative) jumps for 1-shocks (2-shocks) |27j . 

6.1. Nonstationary and Stationary Soliton Edge 

Due to Galilean invariance, there is flexibility in the choice of reference frame for the 
study of DSWs. In the general construction of 1-DSWs presented here, the laboratory 
frame was used with the requirement that the background flow variables lie on the 
1-DSW locus (|5.1a[) . With this coordinate system, three of the four background flow 
properties {pi, p2,ui,U2) can be fixed while the fourth is determined via the 1-DSW 
locus. The soliton and linear wave edge speeds follow according to (|5.10p and (|5.6p 
such that either one of the admissibility criteria for a 1+- or 1^-DSW hold. 

Another convenient coordinate system is one moving with the soliton edge as 
shown in figure [H In this case, one can consider the upstream quantities pi > 0, 
ui > given and the stationary condition 

s±{pi,P2) ^0, (6.7) 

to hold. The downstream density p2 is computed from (j6.7p while the downstream 
velocity U2 follows from the 1-DSW locus ()5.1a|) . The linear wave edge speed is found, 
as usual, from (|5.6|) . 

The admissibility criteria ()6.1a|) - (|6.2c|) for a 1-DSW with stationary soliton edge 
become 

1+-DSW : M2'l<0 < M2 + I, — <Mi- 1, vi < 0, 

Cl 

1--DSW: A/2 - 1< — < A/2 + 1, KMi, U2 > 0. 

C2 

For l+-DSWs, the downstream fiow must be subsonic (A/2 < 1) while for 1~- 
DSWs, the upstream flow must be supersonic (A/i > 1). It is expected that both 
properties hold for both 1^-DSWs but this is not required by the admissibility criteria. 
Supersonic upstream flow and subsonic downstream flow is consistent with classical 
shock waves and will be demonstrated for weak DSWs in section [71 

The DSW locus (|5.1al) , linear wave edge speed in (|5.6p , and the soliton edge speed 
in (|5.10|) along with the stationary condition ()6.7|) constitute the jump conditions 
for a 1-DSW with a stationary soliton edge in dispersive Eulerian fluids. Given the 
upstream Mach number Mi and density pi, the stationary condition ()6.7|) and the 
soliton edge speeds (|5.10p determine the downstream density p2 while the DSW locus 
(|5.1a[) determines the downstream Mach number A/2. 

Figure [1] depicts generic descriptions of 1^-DSWs with stationary soliton edge. 
In the 1~-DSW case of flgure[lja), the soliton edge is the leading edge and the linear 
wave edge is the trailing edge. The opposite orientation is true for the 1+-DSW shown 
in flgure[TJb). This generic behavior was also depicted in [35] based on an analysis of 
weak DSWs in plasma. 
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6.2. Breakdown of Simple Wave Assumption 



A fundamental assumption in this DSW construction is the existence of a simple wave 
or integral curve of the full Whitham modulation equations connecting the trailing 
and leading edge states. This assumption is in addition to the admissibility criteria 
discussed in section |6l The simple wave assumption for a 1-DSW requires a monotonic 
decrease of the associated characteristic speed as the integral curve is traversed from 
left to right |69j . This monotonicity condition leads to the requirement of genuine 
nonlinearity of the full modulation system. Identification of the loss of monotonicity 
is undertaken by examining the behavior of the modulation system at the leading and 
trailing edges. 

The full Whitham modulation system exhibits four characteristic speeds Ai < 
A2 < A3 < A4 that depend on {k,p,u,a). In the case of positive dispersion, the simple 
wave DSW integral curve is associated with the 2-characteristic [55] and connects the 
left, right states (fci, pi, ui, 0), (0, p2, U2, 02), respectively where ki is the wavenumber 
of the wavepacket at the linear wave edge and 02 is the solitary wave amplitude at 
the soliton edge. Generally, the integral curve can be parametrized by p S [pi,p2], 
(fc(p), p, m(p), a(]o)). The monotonicity condition for a simple wave can therefore be 
expressed as 



> 



dA, 



9A, 



k' 



8X2 



(6.9) 



dp dk cTp du da 
where primes denote differentiation with respect to p along the integral curve. At the 
linear wave edge, the A2 and Ai characteristics merge 

z^i(fc,p, m) = lim Ai(A:,p, u, a) = lim A2(fc, p,u,a), (6.10) 

a^-0+ a-)-0+ 

where vi is the smallest characteristic speed of the modulation system when a = 0. 
This merger of characteristics implies that right differentiability of A2 when a = 
requires 



d d 

-X2{k, p,u,0) = — Ai(fc,p,w,0) — 0. 



3.11) 



a oa 

Using (|6.9p . (|6.10p . and (j6.1ip . the breakdown of the monotonicity condition (|6.9p can 
now be identified at the linear wave edge as 



lim 



dA, 



dp 



k' 



ou 



0. 



3.12) 



k=k 



k=kup=pi,u=ui ^ 

The modulation system at the linear wave edge is comprised of the two dispersionlcss 
equations and the conservation of waves |25j 



;i ,p=p\,u=u\ 
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c'/p 
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.13) 



where uj{k, p, u) is the negative branch of the dispersion relation (|2.3p . The eigenvalues 
i^i.2,3 and associated right eigenvectors ri_2.3 for this hyperbolic system are 

(z.i,ri) = (u~ wo., [0,0,1]^), (6.14a) 

(;^2,r2) = (u-c, [p(woj, -c),~c{ujo^ - c),ck~pujOj,f), (6.14&) 

(1^3,^3) = (u + c, [p(wo. +c),c(a;ofc + c),ck -puJoJ'^). (6.14c) 



Shock waves in dispersive Eulerian fluids 



20 



Note the ordering I'l < 1^2 < I's due to the 1+-DSW admissibihty criterion (j6.3&[) . 
Recahing that the 1-DSW integral curve satisfies (|5.3p and (j5.4[) at the hnear wave 
edge, then (|6.12p occurs (breakdown) when 

'^o,*, (^0^+^) +(c-a;oJ (wo,^ + = ) =0- (6.15) 

\ P / V P/ k=ku'p=Pi 

A direct computation shows that dv^/dpi = if and only if (|6.15p holds, offering a 
simple test for linear degeneracy once the DSW speed has been computed. The linear 
degeneracy condition (|6.15p was given in |58| for the Serre equations. Its derivation 
was attributed to loss of genuine nonlinearity in the reduced modulation system ()6.13p 
when a = 0. Linear degeneracy occurs in this system when 

Vl',; • r,; = 0, (6.16) 

for some i G {1,2,3}. 1+-DSW admissibility (j6.36p implies that the 2- and 3- 
characteristic fields do not exhibit linear degeneracy. Evaluation of (|6.16p for the 
1-characteristic field, however, gives 

LOkk{ki,pi) = 0, (6.17) 

corresponding to zero dispersion, a different condition than (|6.15|) . In the vicinity 
of the trailing edge, the 1+-DSW self-similar simple wave corresponds to the first 
characteristic family (|6.14fl|) and satisfies the ODEs 

where differentiation is with respect to the self-similar variable x/t. This demonstrates 
that the Whitham modulation equations exhibit gradient catastrophe, \k'\ — > 00, when 
the dispersion is zero (|6.17p . A direct computation demonstrates that dv+/dp2 = 
if and only if (|6.17p holds. Thus, breaking in the Whitham modulation equations 
coincides with an extremum of the linear wave edge speed with respect to variation in 

Breaking in the Whitham equations has been resolved in specific systems by 
appealing to modulated multiphase waves describing DSW interactions [TOl [71] [72j 
[73l[74]. It is therefore expected that, beyond the zero dispersion point (|6.17|) . internal 
DSW interactions will occur and the simple wave assumption will no longer hold. This 
behavior was intuited by Whitham before the development of DSW theory (see [21] , 
section 15.4) where breaking of the Whitham modulation equations were hypothesized 
to "represent a source of oscillations" . 

An analysis of the soliton wave edge where fc — >■ can be similarly undertaken. 
Recalling that the characteristic speed of the soliton edge is the phase velocity It— (Do/fc 
(|5.10p , the breakdown of the monotonicity condition for the positive dispersion case 
is 

dA2 



lim 

fe-i.o+ dp 



2,P=P2,U= 



dk\k J &p\k 



= 0. 

k — k2 ,p— P2 ^U — U2 



Using the 1-DSW locus (|5.1ap and the characteristic ODE (|5.8p lead to the 
simplification 
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The positivity of the first factor is equivalent to the admissibihty criterion (|6.3ap so it 
is a zero of the second factor 



ck 

^^0-+ — 
P 



= 0, (6.18) 

that offers a new route to hnear degeneracy. Recalhng that the dispersion relation 
involves two branches (|2.3p . care must be taken that the appropriate branch is used 
in (|6.18p . which can change when passing through cjq = 0. A direct computation 
verifies that ds+/dp2 = if and only if (|6.18p holds. Therefore, an easy test for linear 
degeneracy is to find an extremum of s^{pi, P2) with respect to variations in p2- Note 
that the linear degeneracy condition (|6.18p also coincides with a breaking of one of 
the additional sufficient admissibility conditions in (|6.4p . 

Just as zero dispersion at the linear wave edge can lead to singularity formation 
in the Whitham equations, the soliton edge can similarly exhibit catastrophe when 
the phase velocity reaches an extremum 

^ / ~k k=~k2,'P=P2 

This corresponds to zero conjugate dispersion. When (|6.19p is satisfied, wave 
interactions at the leading edge are expected to occur for larger initial jumps. 
In contrast to the linear degeneracy criterion, a direct computation verifies that 
ds+/dpi = if and only if (|5T^ holds. 

The criteria for breakdown of 1^-DSWs is the same as (|6.15p . (|6.17p with 1 — > 2 
and (|6.18p . (|6.19p with 2 1. In summary, two mechanisms at each DSW edge 
for the breakdown of the simple wave assumption have been identified: the loss 
of monotonicity (linear degeneracy) (|6.15p . (|6.18p or gradient catastrophe in the 
Whitham modulation equations due to zero dispersion (|6.17p , (|6.19p . These behaviors 
can be succinctly identified via extrema in the DSW speeds as 

hnear degeneracy — — = or — — = 0, 
1+-DSW: '^P^ <^^2 (g20a) 

gradient catastrophe — — = or — — = 0, 
dp2 dpi 

linear degeneracy — — — or — ^ — = 0, 
1--DSW: dpi (^Q20b) 

dv- 9s_ 
gradient catastrophe — — = or — — = 0. 

dpi dp2 

The negation of these breakdown criteria are further necessary admissibility criteria, 
additional to (|6.1ap - ()6.1c|) . for the validity of the simple wave DSW construction. 

7. Weak DSWs 

The jump conditions for an admissible 1-DSW presented in section [5] apply generally 
to dispersive Eulerian fluids satisfying hypotheses A1-A5. They can be determined 
explicitly in the case of weak DSWs. An asymptotic analysis of the jump conditions 
is presented below assuming a weak 1-DSW corresponding to a small jump in density 

P2 = Pi + A, |A|<1. 
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Two approaches arc taken. First, asymptotics of the Whitham-El simple wave DSW 
closure theory are applied and second, direct KdV asymptotics of the dispersive Euler 
equations are used. 



7.1. Linear Wave Edge 

Expanding the linear wave edge speed (j5.6p yields 

d 

Vj{pi,pi + A) liiR v-j{pup2) + ^v-i{pi,p2)A. (7.1) 

P2-^pi Op2 

Using the long wave asymptotics of the dispersion relation (|2.4p and the initial 
condition for the integral curve (|5.5p . the first term is 

lim vApi,p2) = ui - lim Wq^ = Wi - Ci. 

P2->Pl fc->0 

The derivative term in (17. ip for the case j = 1 is evaluated using (|5.4p - (|5.6p 

lim j—vi{pi,p2)= lim ~uJo^,j— ^ lim wq,,— 
P2^pi op2 P2^pi op2 P2^pi api 

cik/pi + loq- 

The second equality in (|7.2p involves differentiation with respect to the initial "time" 
p2 which, due to uniqueness of solutions to the initial value problem, satisfies 

^iPi;P2) = --^{pi;P2)- (7.3) 
dp2 dpi 

The last equality in (|7.2p follows from the weak dispersion asymptotics (|2.4p . A similar 
computation for the j = 2 case gives 

d dk 
lim j—V2{pi,P2)= Inn u'{p2)-u}o^j,-u}o,,— 
P2->-pi op2 P2^pi dp 

ci A. , ci/c/pi+wo 

= hm ujQ, - + wn, , 

Pi \k^^ " ci-wo. 

Cl 



c 



(7.4) 



Pi 

Combining this result with the other speed calculation gives 

Vjipi,pi+A)^ui~ci + i-iyi3-j) (^j- + c[^A, |A|«1. 

The corresponding wavenumber at the linear wave edge can also be determined 
perturbatively. Note that a Taylor expansion of ki = k{pi] pi + A) for small A is not 
valid because fc(p; pi) is not analytic in a neighborhood of pi, exhibiting a square root 
branch point. However, k^ is analytic, so that upon Taylor expansion, the use of (|7.3p 
and ((^ yield 

3|Ai| \pi 

the wavenumber of the linear wavepacket at a weak 1^'^-DSW's linear wave edge. Note 
that the wavenumber is independent of the sign of dispersion. 



fc.^ — vlAj, J = 1,2, |A|«1, 
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7.2. Soliton Edge 

The soliton edge speed is expanded for a small density jump as 

d 

S]{pi,Pi + ^) ^ lim Sj{pi,p2) + T^Sj{pi,p2)A^ . 

P2^Pl Op2 

Using the expansion p.4|) . the definition (|5.7p . the expression (jS.lOp . and the initial 
condition (j5.9|) gives 



r f \ r i^o{k,Pi) 

hm Sj[pi,p2) = ui — um = = ui — a. 

P2^pi fe_).o k 

To compute the limit -^Sj{pi,pi) necessitates different considerations for each 
j. When j = 1, (|5.10p gives 

9 , X -,. ^Ork- uiQ dk 
hm T^si(pi,p2) ^ lim =r- — 

wof^ - Wo dk 

= hm — — 

P2^Pl /c2 dpi 

(worfc - wo)(cifc/pi + Wo-) 
= hm = 2 

k^a fc2(ci — wo~ ) 

3 \Pi 

When j ~ 2, the limit is similarly computed as 

,. 9 , X r A '^o^ uok-uodk 

hm ;^S2(pi,P2)= hm w (pi) - 1= 

p-2^p\ op2 P2^pi k k^ dp 

- ihn - _ (wo^_A; - wo)(cifc/pi + Wp^) 

fc^o Pi k ^^{cx — Cjq^) 

1 /ci 



3 Vpi ■ 

Combining these results gives the asymptotic soliton edge speed 

s,(pi,pi+A)^Hi-ci-^^ (^^ + c;^ A, |A|«f. (7.5) 

7.5. Admissibility: Positive and Negative Dispersion 

By insertion of the DSW speeds (|7.4p and (|7.5p into the general admissibility criteria, 
it is found that the 1+-DSW criteria (|6.1ap - (|6.1cp are satisfied if and only if A > and 
sgnwfcfe > 0. Similarly, the 1"-DSW criteria (|Oa|) - (p:^ hold if and only if A > 
and sgnwfefe < 0. Hence, the notation 1*-DSW associated with the dispersion sign is 
justified for weak DSWs. 

In the notation of (|6.6p . the weak 1^-DSW speeds are 

4'^(pi,Pi+A)^7.i-ci-^ (^^ + c'i^ A, (7.6a) 
^;£^(pi,Pi + A)^iii-ci-i^ (^^ + c;^ A, 0<A«1, (7.6&) 
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where the superscript denotes the association with a 1-DSW. Notably, the DSW speeds 
(|7.6ap and (j7.66|) differ from the dissipatively regularized shock speed (j4.14p only in 
the numerical coefficient of the (ci/pi + c'i)Sp term. A similar analysis shows that the 
2-DSW locus (j5.16p requires a negative jump in density and yields the speeds 

4'^ (P2 + A, P2) ~ M2 + C2 + ^ ( J + "=2) A, 

V^±\p2+A,p2)^U2+C2 + ^^ (^^+4^ A, 0< A<1. 

7.4- Stationary Soliton Edge 

Choosing the reference frame moving with the 1^-DSW soliton edge so that Sj^-* = 
results in the relations 

A± ^- rT(Mi-l), 

M2,±-l-2±i(il/i-l), 0<il/i-l<l, 

which differ from their classical counterparts (|4.15p only by a numerical coefficient. 
Upstream supersonic flow through a weak, admissible DSW results in downstream 
subsonic flow as in the classical case. 



7.5. KdV DSWs 

An alternative method to derive weak DSW properties is to consider the weakly 
nonlinear behavior of the dispersive Euler equations directly. Inserting the multiple 
scales expansion 

p = pi + ApW (C, T) + A VHe, T) + . . • , 
u = ui- Au(i) , T) + A2u(2) t) + . . . , 
^ = Ai/2[x-(ui-ci)t], T = A3/2i, 0<A<1, 

into (|2.ip . equating like powers of A to ©(A^/^), and recalling the assumed small 
wavenumber behavior of the dispersion relation (12. 4p yields the KdV equation 

41) - (1 + ^) u^'^u^l^ + pug\ ^ 0, = ^u('l (7.7) 

The initial data (|2.5p along the 1-DSW locus (|5.1ap leads to the identification 

z.W(C,0) = J ° ^<'^ (7.8) 
[ ci/pi ^>0. 

The large T behavior of u'^-' satisfying (|7.7p for the initial data (j7.8p results in a DSW 
whose structure and edge speeds depend on sgn p. For negative dispersion, p < 0, the 
DSW is oriented such that the leading, leftmost edge is characterized by a positive, 
bright soliton. The positive dispersion case is equivalent to the negative dispersion 
case by the transformations x — ?> —x, t — > —t, and it*^^-* — > —u'-^\ Therefore, the 
leading, leftmost edge is the linear wave edge and the trailing edge is characterized by a 
negative, dark soliton. Scaling the KdV DSW speeds back to the {x, t) variables results 
precisely in the admissible, approximate weak 1^-DSW speeds (|7.6ap and (j7.6fep . The 
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P2 



x/t 

Figure 3. Universal properties of weak 1^-DSWs with a weak 1-VSW (viscous 
shock wave). DSWs are represented by their envelopes and edge speeds. Speeds 
are in units of (ci/pi + c^)A and ui = ci for simplicity. 

KdV DSW provides additional information, the approximate amplitude of the soliton 
edge 

1+-DSW : p{s+t, t) ^ pi{l~ ^A^ , 
1--DSW : p(s_t, t) ^ pi{l + ^f"^) ■ 

7.6. Discussion 

The analysis of this section has yielded the behavior of weak Eulerian DSWs in 
the context of assumptions A1-A5 and the admissibility criteria (|6.1flP " (|6.2cp . For 
classical, weak 1-shocks (2-shocks), the dispersionless Riemann invariant r2 (ri) is 
constant across the shock to third order in the jump height A [21]. Recalling that 
the DSW loci (|5.1ap . ()5.16p correspond to constancy of a Riemann invariant (simple 
wave condition), the classical Hugoniot loci and the DSW loci for weak shocks are 
the same to 0{A^). However, the shock speeds differ at 0(A). The universal 
properties of weak shocks regularized by positive/negative dispersion and dissipation 
are depicted in figure [3] The jump height A, upstream density pi, and pressure 
law P{p) impart only a uniform scaling of the shock speeds by (ci/pi + c']^)A and 
a relative scaling of the 1^-DSW soliton amplitudes by piA/c^. All shock speeds 
differ, showcasing the distinguishing properties of each regularization type. l~-DSWs 
exhibit backpropagation whereas l+-DSWs and 1-VSWs (viscous shocks) do not. The 
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prominent soliton edge of a 1^-DSW (1+-DSW) propagates faster (slower) than a 
classical shock. By continuity and the discussion of admissibility, it is expected that 
moderate amplitude DSWs for Eulerian fluids with a fixed sign of the dispersion 
exhibit a structure similar to that pictured in figure [1] This is indeed the case for the 
example fluids considered in this work, see section [9] 

8. Dispersive Breaking Time 

In the small dispersion regime, the hydrodynamic dispersionless system (j4.4p 
asymptotically describes the evolution of smooth initial data until breaking occurs. 
One can therefore use the breaking time estimates from section 14.11 to estimate the 
time at which dispersive terms become important. This result was applied to the 
NLS equation with c{p) = p^l"^ in [7S] to estimate the onset of oscillations in fiber 
optic pulse propagation. Here, the breaking time estimates are applied to polytropic 
gases with sound speeds c{p) = p^/^p^/^, < p < 2 (e.g., gNLS with power law 
nonlinearity and gSerre). Generalizations using Lax's theory developed in section |4T1 
are straightforward. 

The flow considered is a slowly varying Gaussian on a quiescent background 

u(x-, 0)=0, /9(a;, 0) = 1 + exp[-(ea;)^], < £ < 1. 

Slowly varying initial data ensures the applicability of the dispersionless system (14. 4p 
up to breaking when t ~ OiXje). This choice of initial data has been used in photonic 
DSW experiments [33]. 

Figure|3]shows the results for the gNLS equation (j3.ip with power law nonlinearity 

/(p) = ■ The solid curves ( ) correspond to the upper and lower bounds on 

the dispersionless breaking time estimates ()4.11|) and (• ) correspond to numerically 
computed breaking times for several choices of p and e. The breaking time from 
simulations is defined to be the time at which the density first develops one full 
oscillation in the breaking region. For \e\ <SC 1, the bounds (|4.1ip accurately estimate 
the breaking time across a range of nonlinearities p. Note that, for these dispersive 
Euler models and class of initial data, e ^ 0.01 is required to obtain agreement with 
the dispersionless estimates. 

9. Large Amplitude gNLS DSWs 

The general Whitham-El simple wave DSW theory is now implemented for the gNLS 
equation. 

To simplify the presentation, and without loss of generality, the independent and 
dependent variables in (|2.2p will be scaled so that the initial jump in density (|2.5p is 
positive, from unit density 

Pi = 1, P2 = A > 1, 

so that 1-DSWs will be considered. Then, according to the 1-DSW locus (|5.1ap of 
admissible states, the jump in velocity satisfies 
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Figure 4. Breaking time bounds ( ) and numerically computed breaking 

times (• ) for gNLS with power law nonlinearity /(p) = p^. 



9.1. General Properties 

The gNLS equation exhibits positive dispersion so that DSWs are of the l"*" variety. 
Substituting the expressions p.8p and (|5.3p into the ODEs (|5.2p and (|5.8|) resuhs in 
the initial value problems 

1 



1 

dfc p 



4p/'(p) 



1 

2p 



f"{-p) 
2/'(p) 



1 



1 



1/2 



fc(A)=0. (9.1) 



2p/'(p) L 4p/'(p) 
for the determination of the linear wave edge speed and 

1/2 



dfc 
d^ 



1 

-fc — 



1 



4p/'(p) 



1 
2p 



2/'(p) 



1 



2p/'(p) 



1 



4p/'(p) 



1/2 



fc(l) - 0, (9.2) 



for the soliton edge speed. In (|9.2p . it is possible for the quantity within the square 
roots to pass through zero. When this occurs, an appropriate branch of the dispersion 
relation should be used so that the conjugate wavenumber remains real valued. 
Recalling remark [T] in section 15.21 the transformation 

1/2 



a{p) 



uJo{k,p) 



c{p)k 

simphfies dH]) to the ODE 

da 1, 
d^^-2(l + ' 



1 



4p/'(p) 

(2a - l)/"(p) 
(2a+l)/'(p) 



with initial condition 
a(A) 



1. 



(9.3) 



(9.4) 
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The analogous transformation for the conjugate variables 



a(p) 



Wo(A:,p) 



1 



1/2 



transforms 



4c(p)2 

with a ^ a and the initial condition 



5) 



c{p)k 

to the same equation 

a(l) = 1. (9.6) 

Upon solving the initial value problems for a and a, the linear wave and soliton 
edge speeds are found from (|5.6p and (|5.10p . respectively, which take the form 

l-2a(l)2 



u+ = Ui + 



a{l) 



c 






. p . 



.1/2 



dp, 



(9.7a) 



(9.7&) 



in the a, a variables. The weak DSW results (|7.4p and (|7.5p give the approximations 



v+ ^ Ui 



vm 



1 



ni) 
/'(I) 
, /"(I) 



/'(I) 



(A-i; 



(A-l) 



< A- K 1. 



(9.8a) 
(9.86) 



Equating the soliton speed in (j9.76P to the soliton/amplitude speed relation in 
(|3.10p gives an implicit relation for the dark soliton minimum p,nin in terms of the 
background density A 

r-A 



.(A)2A/'(A) = 



2pn 



Pmin 



/(A)- 



1 



A 



Pmin Jp, 



./(P)dp 



(9.9) 



A direct computation shows that neither linear degeneracy (j6.15p nor a sign of 
dispersion change (|6.17p occurs at the linear wave edge. However, at the soliton 
edge, there are several distinguished values of a(A) with physical ramifications. The 
modulation theory breaks down due to singular derivative formation in (j9.3p when 

a(A = A,) = -i. 

From the denominator in (|5.8p . singularity formation occurs precisely when cDq^ = 
c(A). A direct computation shows that luq is a concave function of fc, which implies 
ujo/k > = c(A) so that singularity formation coincides with the violation of the 
admissibility criterion (|6.3ap . 

From the initial condition (|9.6p and the ODE (|9.3p . a{A) decreases from 1 for 
increasing A sufficiently close to 1. The value of a at which its derivative is zero, from 
([931) . is 

A/"(A) - /'(A) 



.(A) 



2[A/"(A) + /'(A)]- 



So long as 



a{A) > max 



.(A),-i 



the right hand side of (19.31) is strictly negative and finite. Since a satisfies the same 
ODE as a, if ciniin < 1, ck(1) is an increasing function for all A and q;(A) is a decreasing 
function of A until a{A) = max(Q;niin(A), —1/2). 
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As A is increased from 1, cavitation (a point of zero density or vacuum point) first 
occurs when Pmin ~ corresponding to a black soliton at the traihng edge, moving 
with the background flow speed u(A). Since necessarily A > pmin, (|9.9|) implies that 
a vacuum point first occurs when 

a(A = Av) = 0. (9.10) 

For larger jumps A > Av, the vacuum point is expected to develop in the interior of 
the DSW [Ml ED]. According to the transformation (|9.5|) . when a crosses zero, the 
branch of the conjugate dispersion changes sign. 

Using Galilean invariance, it is convenient to consider the reference frame moving 
with the soliton edge so that s+ = 0. In such a frame, given the upstream supersonic 
flow velocity ui > 1, the density jump A is determined from (|9.76P and the downstream 
flow velocity satisfies U2 = a{A)y/ Af'{A). Then, for A < Av, "(A) > so that 
U2 > 0. But, when A > Av, cavitation occurs and (5(A) changes sign causing U2 < 0. 
Counterintuitively, the dispersive fluid flows into the DSW from both sides upon the 
generation of a vacuum point. This behavior has been observed in NLS [24] and 
photorefractive gNLS [50] . 

At the soliton edge, linear degeneracy (|6.18p occurs when 



ck ~ 

LUp- + — = ck 

P 



f'{A) + Af"iA) ^ 1 



= 0. (9.11) 



2A/'(A)<5 A_ 

The only way for this to occur is for ujq to change to another branch of the conjugate 
dispersion relation, i.e., for a to pass through 0. Then from (|9.1ip . the value of 5 at 
which linear degeneracy occurs is 



5i(A) = -i 



^ , A/"(A) 



(9.12) 



/'(A) 

According to the assumptions on / p.7|) . fii < as required. Since the linear 
degeneracy condition (|9.1ip amounts to ds+/dA = (recall (|6.20a|) ). the distinguished 
value ai coincides with an extremum of the soliton edge speed as the jump height A 
is varied. 

A direct computation verifies that the zero dispersion criterion (|6.19p does not 
occur. The admissibility criteria (|6.1aP and (|6.16p correspond to a{A) < 1 and 
a(l) < 1 which are true generically. In summary, so long as cimin < —1/2, (5(A) 
is a decreasing function that can attain the distinguished values (5 = when A = Av 
corresponding to a vacuum point or cavitation, a = ai when A = Ai coinciding with 
linear degeneracy and the breakdown of the simple wave criterion, and a = —1/2 
when A = As leading to singularity formation. For A < Ai, the only admissibility 
criterion left to verify is the leading and trailing edge ordering (|6.1cp . 

The DSW regularization is completed upon solving the ODE (|9.3p . This equation 
is separable for nonlinearity satisfying f'{p) ± pf"{p) = 0. Recalling the admissible 
nonlinearity p.7p . the cases of interest are f{p) = npP , p > 0, corresponding to 
polytropic superfluids. This class of nonlinearity is now considered. 

9.2. Power Law Nonlinearity 

With the choice p.2p of power law nonlinearity, the speed of sound is c{p) ~ p^^^pP^^. 



Equation (|9.3|) becomes the separable equation 



-2(1 + 2") -d«^?. (9.13) 



{l + a)i2-p + 2pa) p 
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Using the initial condition (j9.4[) . (|9.13p is integrated to 



l + a(l) 



2-p + 2pa{l) 



2(p-l)/p 



^ ^(3p-2)/2^ 



2 In 



+ 



1 



l + a(l) 



P = n 



and for ()9.6p . (|9.13p is integrated to 
l + a(A) r2-p + 2p(5(A) 



2(p-l)/p 



^(2-3p)/2 



2 In 



l + a(A) 



+ 



1 



l + a(A) 2 3 



(9.14a) 
(9.14&) 

(9.15a) 
(9.15&) 



the integrals existing so long as 



p~ 2 

a(l), a(A) > — — = a^nin- 

The linear wave and soliton edge speeds are then determined via (|9.7ap and (|9.76p 
w+ - ui _ 1 - 2q;(1)2 



Cl 

s+ - Ui 

Cl 



«(1) 



= -AP/2 



.(A) 



Equations (|9.14ap and (|9.15aP can be solved explicitly for a in several cases 



25^1/4 3 



'^(1) 16- 2 

.(A) = -A-i/4 - - + A 
^ ^ 16 2 16 



25A1/2 _ leAi/-* 



1/2 



25A-i/2_i6A-i/4 



1/2 



p=l 



p = 2 



a(l) = 2Ai/2_i, 
5(A) = 2A-1/2 _ 1^ 

a{l) = \ fv/l + 8A2-l) , 



a(A) = i(Vl + 8A-2-l), 



(9.16a) 
(9.166) 

(9.17a) 
(9.176) 
(9.17c) 



providing explicit expressions for the 1-DSW linear wave and soliton edge speeds 

2 



p=l 



Cl 

S+ - Ul 



16-i 



25A1/4 _ 24 + 5V25AV2 _ i6Ai/4 



Cl 
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Cl 

s+ - 


Ul 


Cl 






Ul 
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Ul 



_5^i/4 39 _ 1^25-16Ai/4, 
2 16 16 



25A1/4 - 24 + 5V25A1/2 _ l6Ai/4 

1/4 I ^ [ 

16 1 

8A - 8A1/2 + 1 
2A1/2 _ 1 ' 

-Ai/2, 



= -2 



Cl 



4A2 - VI + 8A2 
\/l + 8A2 

i (a - 2 + V8 + A2) 



(9.18a) 



(9.186) 



(9.18c) 
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The cubic NLS case (|9.186p agrees with the original result in [24] . 

The distinguished values of the jump height A predict different DSW behavior 
as p varies. The large A behavior of (I9.15al) and the fact that a (A) is a decreasing 
function of A proves that 

ci(A)\^>-i A^oo, p>l, (9.19) 

2p 2 

precluding the possibility of singularity formation when p > 1. Furthermore, when 
p> 2, (|9.19p prohibits cavitation in the DSW because a(A) > 0. 

The cavitation condition (|9.10p can be used along with (j9.15ap to determine 
Av(p), the value of A at which a vacuum point is predicted to appear 

A.(p)= 2^/(3^-2) (^l^j , 0<p<2, p^l, (9.20a) 

For p approaching 2, the vacuum jump Av increases without bound. The limiting 
behavior for p ^ is A^{p) e'^/2 sa 3.69. 

The DSW soliton speed exhibits a minimum (maximum in absolute value) for 
< p < 1 corresponding to the onset of linear degeneracy in the modulation system. 
The value of the jump height Ai at the linear degeneracy point is found by inserting 
([O^ into (|9.15a|) and solving for A 

2/(2-3p) 

Mp) =" 



1(2 - p)ii - pr^/^ 



0<P<1, (9.21a) 
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The p — > behavior is limp_>o Ai(p) = e^/2, coincidentally the same limiting value 
Av(0). The maximum value of the trailing edge speed |s+| is therefore 



max s+ = ui + ci 

l<A<As 



0<p<l. 



2p p 

For A > Ai(p), the simple wave construction is no longer valid. 

For < p < 1, 5 (A) exhibits singularity formation in its derivative when 
a(A) -1/2 for A As{p) where 

A,(p) = 24/b(3p-2)]L^ ^ 0<p<l, p^2/3, (9.22a) 



As(2/3) = 64exp(-9/4) w 6.75, (9.226) 

so that DSWs with A > As{p) are not admissible. 

Based on the foregoing analysis and upon verification of the DSW edge ordering 
(|6.1cp . the modulation theory presented here does not violate any of the admissibility 
criteria for all 1 < A when p > 1 and for 1 < A < Ai(p) when < p < 1. 

Utilizing the special jump values Ag, Av, and Ai, figure [5] right shows regions 
in (p. A) space where a vacuum point is predicted to develop, the DSW modulation 
theory breaks down due to linear degeneracy, and singular behavior. It is clear from 
this figure that Av < Aj < Ag so that the simple wave condition has been broken due 
to linear degeneracy before singularity formation when A = Ag. 
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Figure 5. Right: singular A > As(p), linearly degenerate As(p) > A > Ai(p), 
and cavitation Ai(p) > A > Av(p) regions for a simple wave led gNLS 1+-DSW 
with power law p nonlinearity. Left: numerically computed l+-DSWs for p = 2/3. 



The simple wave DSW construction makes several predictions that can be tested 
through long-time numerical simulations of the shock tube problem of section 15.31 
For the details of the numerical method, see the appendix. Figure 6(a) depicts the 
numerically computed density at t ~ 350 resulting from an initial jump in density 
(pi,P2) = (1,14) at a; = 2000 for power law gNLS with p = 2/3. The structure of 
a 1+-DSW connected to a 2-rarefaction via an intermediate, constant state is clear. 
Almost imperceptible modulations of the zoomed-in solution in the inset demonstrate 
the two scale nature of dispersive hydrodynamics. Figure [6 (b ) | provides a comparison 

of the intermediate density Pm predicted by (j5.13p ( ), numerical computation (■, 

A, • ), and the intermediate density predicted for a classical, polytropic gas, equations 

(|4.17p ( ). The integrable NLS case p ~ 1 agrees precisely with the numerical 

computations but for the non-integrable cases, eventual deviation is observed for large 
initial jumps. Interestingly, the classical shock prediction begins to agree with the 
p = 2 numerics for sufficiently large jumps. Similar behavior was noted for gNLS with 
photorefractive nonlinearity |50j . 

Returning to figure [5] left, numerical simulations for p = 2/2) and varying A are 
displayed. A transition to cavitation is shown for sufficiently large jumps. While the 
DSW modulation theory breaks down due to linear degeneracy for a larger jump, the 
simulations in figure [S] show no noticeable change in the computed DSW structure. 

A plot of the resulting 1+-DSW speeds for several values of p, normalized by the 
sound speed ci , is shown in figure [T] The value of A corresponding to the numerical 
simulation results (■, A, • ) is the one extracted from the computation. The A value 

used for the theoretical predictions ( , ) is pm from the shock tube problem 

(|5.13p . The computations for the integrable case p = 1 agree excellently, even for large 

jumps. The solid curves ( ) show rapid deviation from the weak DSW straight 

line ( ) predictions, signaling the important effects of nonlinearity. In the weak to 

moderate jump regime A < 3, the non-integrable cases exhibit good agreement with 
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Figure 6. (a) Density resulting from the numerical solution of the shock tube 
problem for power law gNLS. The inset depicts the slowly varying nature of the 

DSW. (b) Predicted ( ) and computed (■, A, • ) intermediate states pm 

resulting from the shock tube problem. The prediction for a classical shock tube 
is also given ( ). 



the predicted soliton edge speed but deviation occurs for large jumps. Recall that for 
p = 2/3, linear degeneracy sets in for jump heights above that which corresponds to the 
minimum of the soliton edge speed curve (|9.215p . Ai w 6. However, the computed 1'''- 
DSW shows no distinct change in structure for jumps below and above this threshold. 



both looking qualitatively similar to the DSW depicted in figure 6(a) Interestingly, 



the non-intcgrable cases admit excellent agreement with the predicted linear wave 
edge speed across their regions of validity. 




Figure 7. Power law gNLS 1-DSW speeds for varying density jump A and 
nonlinearity exponent p. Solid curves ( ) arc from ||9. 1661 1. 119. 16 al l. (■, A, • ) 

result from numerical simulation, and ( ) arc the weak DSW results l|9.8a|l . 



The soliton minimum pmin, determined by solving (|9.9p with a numerical 
rootfinder, is shown in figure [5] for different choices of the nonlinearity p and jump 
heights (curves) . The simple wave DSW construction predicts no cavitation for p > 1 
and a sharply increasing soliton minimum as the jump height is increased beyond 
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Figure 8. Power law gNLS DSW trailing soliton minimum Pmin as a function of 
the DSW jump and p. 



Av(2/3). Neither of these behaviors are observed numerically. Sufficiently large jumps 
in the p = 2 case lead to cavitation. The numerically computed density minima of the 
trailing edge soliton (■, • ), curiously, lie on the integrable NLS curve. It is not 
clear as to why this is the case. 



9.3. Nonpolynomial Nonlinearity 

Another particular form of the nonlinearity is now considered, that of the 
nonpolynomial nonlinearity p.3p . When inserted into (|9.3[) . this nonlinearity results 
in the non-separable ODE parametrized by 7 

dg^ (l + a)[2 + 37p + 2a(2 + 7p)] 

dp 4p(l+7p)(l + 2a) ■ ^ ' 

It is necessary to numerically solve (|9.23p with the initial conditions (|9.4p . or (|9.6p 
with a — >■ a, in order to recover DSW properties. Since 

a.ni„(A) = -^ + ^<-i 7A>0, 

and a(A) is a decreasing function for a > amin, singularity formation (a — )- —1/2) 
is guaranteed for 7 > and A sufficiently large. Linear degeneracy, however, occurs 
before this when 

a(Ai) = — > — . 

^ ^' 4 + 47A1 2 

Note that for 7^1, the characteristic equation (|9.23p asymptotes to the power law 
nonlinearity equation (j9.13p with p ~ 1/2 so that the behavior in the large 7 limit can 
be recovered directly from (|9.17a|) and (|9.18ap . 

Similar to the power law nonlinearity, the three distinguished values for the shock 
jump Av < Ai < As divide (7, A) parameter space into various regions as shown in 
figure [HI The large 7 behavior of the distinguished values are, from (|9.20ap , (|9.21ap , 
and (|9.22ap with p=l/2 

Av(7)~3.7, Ai(7)^5.1, As(7) ~ 6.0, 7 » 1- 

Cavitation and linear degeneracy are predicted to occur for 7 > with sufficiently large 
A. Additionally, the small 7 asymptotics can be determined as a regular perturbation 
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Figure 9. Singular A > As(7), linearly degenerate As(7) > A > Ai(7) and 
cavitation Ai(7) > A > Av(7) regions for a simple wave led 1+-DSW with 
nonpolynomial nonlincarity. 
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Figure 10. Predicted ( ) and computed (■, A, • ) intermediate states pm 

resulting from the nonpolynomial gNLS shock tube problem. 



to the NLS solution ()9.17&|1 . The result is complex so it is not reported here, deferring 
rather to the numerical solution of the ODE (|9.23p . 

The shock tube problem for 7 € {0, 1/4, 1} was solved numerically. The resulting 
DSW shock structure exhibits a qualitatively similar pattern to the one shown in 



figure 6(a) The computed intermediate density prn as compared with that resulting 
from the 1-DSW locus (|5.13p is shown in figure [101 As before, when deviating from 
the integrable case with 7 > 0, the 1-DSW locus provides a good prediction to the 
intermediate density for weak to moderate jumps but deviates in the large jump 
regime. 

The computed DSW speeds arc pictured in figure [11] for varying 7 and A. The 

weak DSW results ( ) rapidly deviate from the predicted ( ) and computed 

(■, ▲, • ) results for A > 1.5. The non-integrable DSW speeds with finite 7 deviate 
from the DSW regularization in the large jump regime A > 3, more so for the soliton 
edge speeds. The onset of linear degeneracy and singular behavior is clear from the 
solid soliton speed curves, exhibiting minima at A = Ai(7) and an infinite derivative 
at A = A,(7). 

The soliton minimum pmin is shown in figure 1121 As with the power law 
nonlincarity, the computed soliton trailing edge minima lie approximately on the 
integrable NLS curve suggesting that this behavior is not coincidental to a particular 
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(a) Soliton edge speeds (b) Linear wave edge speeds 

Figure 11. Nonpolynomial gNLS 1-DSW speeds for varying density jump A 
and parameter 7. Solid curves are from l|9.166|l . H9.16a|l . {■, A, • ) result from 




Figure 12. Soliton minimum Pmin as a function of the DSW jump and 
nonpolynomial nonlinearity coefficient 7. 

type of nonlinearity chosen. 
9.4- Photorefractive Media 

The case of DSWs in photorefractive media p.4p was studied in detail in [5D]. There, 
the modulation theory was found to diverge for sufficiently large amplitude jumps and 
an explanation in terms of the existence of cavitation was given. Based on the theory 
developed in section 16. 2[ this breakdown can be identified with linear degeneracy 
at the soliton edge, hence the simple wave condition no longer holds. A tell-tale 
sign of breakdown is realized by the extremum of the shock speed with respect to 
jump height in figure 7 of [5D]. Figure [T3] depicts the phase space (7, A) divided 
into regions by the curves Ai(7), Av(7), and As(7). In contrast to the previous 
discussions, photorefractive nonlinearity only preserves the ordering Av < Ai < Ag 
for 7 ^ 0.287. However, the generic nonlinearity assumptions p.7p corresponding to a 
convex pressure law or, cquivalently, an increasing sound speed are violated for jumps 
exceeding A > I/7. The dashed curve A = I/7 shown in figure [T51 demonstrates that 
the change in ordering of Ai and Av occurs when c'(A) ~ 0. 
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Figure 13. Validity of DSW modulation theory for photorefractive nonlinearity. 
The dashed curve corresponds to c'(A) = 0. 

10. Large Amplitude DSWs with Negative Dispersion 

Dispersive shock waves for fully nonlinear shallow water waves and ion-acoustic plasma 
have been studied elsewhere |25[ I58j . It is worth commenting on these results in light 
of the theory presented here. Both systems exhibit alternative negative dispersion 
regularizations of the Euler equations in contrast to the positive dispersion of gNLS. 

In shallow water with zero surface tension a = |58| . the linear degeneracy 
condition (|6.15p occurs at the trailing edge for relatively small jumps A = Ai « 1.43. 
For A > Ai, the simple wave theory breaks down and, indeed, the numerically 
computed soliton edge speed begins to noticeably deviate from its theoretical value 
(see figures 3 and 4 in [58]). It is curious that the numerically computed linear edge 
speed remains fairly accurate. For a = 0, zero dispersion does not occur at either 
edge, nor does linear degeneracy at the soliton edge. For cr > 0, both uJkk = and 
{Cj/k)j, = when p = \/3a, offering a potential route to gradient catastrophe in the 
Whitham modulation equations for appropriate DSW jumps. 

Ion-acoustic plasma DSWs also exhibit linear degeneracy at the linear wave edge 
for appropriate jump heights. The linear edge speed plotted in figure 5 of [55] for 
a 2~-DSW from pi = A to p2 = 1 exhibits the requisite minimum as A is varied 
(recall (|6.20&[) ). A numerical computation shows that the minimum occurs when 
A = Ai « 1.41. Zero dispersion at either edge and linear degeneracy at the soliton 
edge do not occur by calculation of the criteria from section 16.21 

The particular dispersive Eulerian fiuids discussed here demonstrate the 
importance of the breakdown criteria (j6.20ap . (|6.20&p in both negative and positive 
dispersion cases. 

11. Conclusion 

Shock waves in the dispersively regularized isentropic Euler P-system under modest 
assumptions were constructed using the Whitham-El simple wave closure technique. 
A complete, explicit characterization of admissible, weak DSWs was shown to depend 
only on the pressure law and convexity or concavity of the dispersion relation. Linear 
degeneracy of the modulation equations and zero dispersion leading to gradient 
catastrophe in the modulation equations were identified as causes of the breakdown 
of the simple wave assumption. Simple tests in terms of extrema of the DSW leading 
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or trailing edge speeds for these behaviors were elucidated. Large amplitude DSWs 
were constructed for the case of the gNLS equation modeling super and optical fluids. 
Comparisons with careful numerical simulations of the shock tube problem reveal 
excellent agreement with theory in the weak to moderate jump regime. Deviation 
occurred in the large jump regime with linear degeneracy in the modulation equations 
at the DSW soliton edge a proximate cause for certain pressure laws. 

While there are a number of parallels between classical, viscous Eulerian fluids 
and dispersive Eulerian fluids, this work has demonstrated that DSWs exhibit 
distinct physical and mathematical behavior. Physically, the generation of oscillations 
leads to an expanding oscillatory region with two speeds in contrast to localized, 
classical shock fronts propagating as traveling waves. Positive dispersion fluids 
can exhibit DSWs with cavitation while negative dispersion fluids admit DSWs 
with backflow. Mathematically, weak Eulerian DSWs and shocks exhibit universal 
behavior, depending only upon the sign of dispersion and the pressure law. In 
the large amplitude regime, universality is maintained for classical viscous shock 
jump conditions, which are the same for a large class of dissipative regularizations. 
However, due to their nonlocal nature and the weak limit involved, large amplitude 
DSWs crucially depend upon the particular form of the dispersion. Furthermore, 
admissibility of simple wave led DSWs is much more subtle than the elegantly stated 
Lax entropy conditions for classical shocks. 
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Appendix: Numerical Methods 

The mmrerical solution of the gNLS equation for the shock tube problem, the 

initial step in density 



/P2 X> Xq 

is briefiy described here. A pseudospectral, time-splitting method is implemented for 
the accurate solution of long time evolution for x G (0,i). The initial data (|A.1[) is 
smoothed by use of the hyperbolic tangent initial condition 

r 1 ^ 

V'(a;,0) = I -[l + tanh(a:;o-x)](l-p2) + /02 

where xq = L/2. Time stepping proceeds by use of second order Strang splitting |76j 
where the linear PDF 

is advanced half a time step At/2 exactly followed by a full time step of the nonlinear 
ODE 



dt 



/(IV'NLr)V'NL, VNL(2;,i) = V'L(a;,t + At/2). (A.3) 
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The linear PDE is then advanced half a time stop with the initial data ipLi^, t+At/2) = 
tpNhiXjt + At) giving the second order accurate approximation of ip{x,t + At) w 
ipi^{x,At). Equation (|A.2[) is projected onto a truncated cosine basis of N terms 
that maintains Neumann (ip^ = 0) boundary conditions, computed efficiently via the 
FFT, and integrated explicitly in time. The nonlinear ODE (jA.3[) conserves IV'nlP 
so is also integrated explicitly in time. The parameter Aa; — L/N is the spatial 
grid spacing of the grid points xj = Ax{j — 1/2), j = 1, 2, . . . , A^. The accuracy of 
the solution is monitored by computing the relative deviation in the conserved 
norm E{t) = {ipix, t)\'^dx, E.^i = \E(tf) - E{0)\/E{0) where tf is the final time. AU 
computations presented here exhibit E^^i < 10~®. Also, the accurate spatial resolution 
of the oscillatory structures is guaranteed by the fact that the coefficient of the largest 
wavenumber in the cosine series is less than 5 • 10~^° (often times much less). The 
numerical parameters L, At, and t[ vary depending upon the nonlinearity strength 
and jump height. For example, for power law gNLS with p = 2 and p2 > 11, A^ = 2^^, 
L ^ 1200, At = 0.0002, and ti = 30. Whereas, for p = 2/3 with p2 = 2, N ^ 3 ■ 2^^, 
L = 3000, At = 0.002, and ti = 500. 

The extraction of the DSW speeds s+, and minimum density pmin is performed 
as follows. The precise location of the DSW soliton trailing edge is computed by 
creating a local cubic spline interpolant through the computed grid points in the 
neighborhood of the dark soliton minimum. A root finder is applied to the derivative 
of this interpolant in order to extract the off-grid location of the soliton edge Xs{t) 
and /Omin = |V'(2^s(if), if)P- The slope of a linear least square fit through Xs{tj) for 
j = 1,...,100 equispaced tj G [tf — l,tt] determines s+. For the leading, linear 
wave edge, an envelope function is determined by least squares fitting two lines, each 
through about 30 local maxima and minima, respectively, of the DSW density in the 
vicinity of the trailing edge. The extrema are computed the same as for the soliton 
minimum. The point of intersection of these two lines is the location of the linear 
wave edge Xi,{t). The same fitting procedure as was used to determine s+ from Xs{t) 
is used to extract w+ from Xjj{t). 
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